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We simulate spectral functions for electron-phonon coupling in a filled band system - far from the 
asymptotic limit often assumed where the phonon energy is very small compared to the Fermi energy 
in a parabolic band and the Migdal theorem predicting (1 + A) quasiparticle renormalizations is valid. 
These spectral functions are examined over a wide range of parameter space through techniques 
often used in angle-resolved photoemission spectroscopy (ARPES). Analyzing over 1200 simulations 
we consider variations of the microscopic coupling strength, phonon energy and dimensionality for 
two models: a momentum-independent Holstein model, and momentum-dependent coupling to a 
breathing mode phonon. In this limit we find that any 'effective coupling', Aofi, inferred from the 
quasiparticle renormalizations differs from the microscopic dimensionless coupling characterizing 
these Hamiltonians, A, and could drastically either over- or under-estimate it depending on the 
particular parameters and model. In contrast, we show that perturbation theory retains good 
predictive power for low coupling and small momenta, and that the momentum-dependence of the 
self-energy can be revealed via the relationship between velocity renormalization and quasiparticle 
strength. Additionally we find that (although not strictly valid) it is often possible to infer the self- 
energy and bare electronic structure through a self-consistent Kramers-Kronig bare-band fitting; 
and also that through lineshape alone, when Lorentzian, it is possible to reliably extract the shape 
of the imaginary part of a momentum-dependent self-energy without reference to the bare-band. 



PACS numbers: 71.38.-k, 79.60.-i, 74.25. Jb 
INTRODUCTION 

The many-body problem allows relatively simple in- 
teractions to transform into a wide range of exciting yet 
often complicated phenomena. The quasiparticle pic- 
ture simplifies these complications by grouping funda- 
mental particles and excitations together into quasipar- 
ticles, which themselves behave in a more understandable 
manner. In this picture the real part of the self-energy 
represents the energy difference from the bare particle 
energy, and the imaginary part the inverse lifetime of the 
combined excitation. Angle-resolved photoemission spec- 
troscopy (ARPES) is a well established tool for the inves- 
tigation of such electronic systems as it provides access 
to the electron-removal part of the momentum-resolved 
spectral function A(k, w) [T], which is generally written 
in the form: 

, , _1 S"(k,cj) 

^ 7r[cj-e^-I]'(k,w)]2 + [I]"(k,c^)]2- ^> 

The analysis of this extremely rich data source can be 
both difficult and rewarding as it depends on both the 
interaction self-energy E(k, w) = I]'(k, cj) -|-iE"(k, w), as 
well as the single-particle electronic dispersion (the so- 
called 'bare-band'). A variety of approaches to analyzing 
this spectral function have been utilized, and often focus 
on analysis of either quasiparticle dispersions and their 
path through (k, oj) space, or lineshape and its implica- 
tions for the structure of the self-energy. Both methods 
generally cut the spectral function into curves constant 
in either momentum (generating a series of energy distri- 
bution curves [EDCs]) or energy (for a series of momen- 



tum distribution curves [MDCs]). In this work, using 
simulations which have no experimental limitations, we 
will perform quasiparticle analysis on EDCs (which al- 
lows the identification of a quasiparticle peak in each k 
slice thereby forming a quasiparticle dispersion, e^) and 
self-energy analysis on MDCs (as self-energies often show 
stronger energy dependence, allowing the approximation 
of a constant value over a slice of constant energy). 

In quasiparticle analysis one can estimate properties 
such as the dispersion's velocity = 9£^/9k, eflFective 
mass m^, where 1/m^ = d^e^/dlc^, and quasiparticle 
strength Z^, where Zl = pA{'k,uj)d uj is the integral over 
the coherent part of the spectral function (this is the 
quasiparticle weight only, which in a somewhat loose ter- 
minology is often referred to as quasiparticle coherence 
Ij). If the bare-band dispersion is known, the renor- 
malization of these properties can also be calculated. 
This concept has been used to generate an 'effective cou- 
pling' (which we will denote Acff, but which is often de- 
noted simply A in ARPES literature) in the analysis of 
many complex systems, often through the so-called 'mass 
enhancement factor' m^/m^ = v^/v^ ^ ^ l/(l+Aoff). 
This factor has become a de facto standard in ARPES 
analysis [IH5] since, in the Migdal/Eliashberg limit after 
few approximations, it is equivalent to the true dimen- 
sionless microscopic coupling found in the Hamiltonian 
(denoted A here and in theoretical literature) and is ex- 
pected to manifest itself in a variety of different mea- 
surements |Tj. However, the large values sometimes 
measured for these renormalizations and effective cou- 
plings via ARPES (see, for example, Refs.lMTT|). should 
raise the question of this scheme's universal utility \12\ - 
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FIG. 1: (Color online). The spectral function (a-c) and self- 
energies (d-f) for the ID momentum-independent Holstein 
polaron model calculated with MA'^-* for f2 = 50meV and 
different microscopic couplings, with the bare-band e| , quasi- 
particle band el, and fcm(<^) path shown. 

116] . and generally whether the limits implied by such 
analysis do apply to all systems being measured [T7ti20j . 

Another common goal of spectral function analysis is 
to extract the self-energy. In most circumstances, un- 
der the assumption of k-independence of the self-energy, 
MDCs cuts through Eq.[l] reduce to a simple Lorentzian 
form, thus allowing a measurement of S'(aj) and Yi"{uj) 
through ARPES O [El EIHIS]. However not only do 
these methods hinge on some assumptions and/or ap- 
proximations for the bare-band e^, but more funda- 
mentally the problem of how momentum-dependence in 
S'(k, cj) and E"(k, w) affects this analysis is unaddressed 
- even though it is known that a Lorentzian line shape 
does not guarantee a momentum-independent self-energy 

m- 

Here we present a methodological study of established 
methods and present some new variations using one of 
the most studied interactions - that of electrons and 
phonons. We generate self-energies and spectral func- 
tions where the inclusion of momentum-dependence and 
all energy scales are controlled using the least compli- 
cated electron-phonon interaction models possible. How- 
ever, these models lie outside the limits of Migdal's the- 
orem ^7] where the Eliashberg textbook definition of 
(1-l-A) renormalization is expected to be valid. Before we 
delve into our findings for quasiparticle and self-energy 
analysis, we will first illustrate some aspects of our chosen 
models and how they are simulated. 



THE MODELS 

We use single electron addition to an empty band 
to simulate photoemission from a completely filled sys- 



tem, at Kelvin. Note that this case can be exactly 
mapped onto that of a single particle added to an empty 
band through particle-hole symmetry, which essentially 
amounts to replacing w — > — w. This is an ideal test 
case as it provides the simplest possible description of 
electron-phonon coupling and is uncomplicated by fur- 
ther interactions such as strong correlations between elec- 
trons (as in, for example, Ref.f28l). or even a Fermi sea 
which would add yet another energy scale to the problem 
(as in, for example, Refs. l^ I5D|) . The chemical potential 
in our treatment is then the top of the first electron re- 
moval state, labelled as binding energy on all plots. For 
momentum-independent study we will use spectral func- 
tions and self-energies generated with the momentum- 
average approximation MA^^^ [31] [35]. Since MA(i) has 
been shown to be accurate everywhere in parameter space 
[32] . it will enable us to study A(k,uj) and S(a;) over a 
broad range of electron-phonon coupling and phonon en- 
ergies. For momentum-dependent study we use an ex- 
tension of this approximation with variational considera- 
tions, denoted MA^"'"^ [33]. Although generally accurate 
everywhere in parameter space, for reasons specific to 
this approximation details studied through EDC quasi- 
particle analysis are best realized through MA^'"'^\ and 
MDC based self-energy analysis is best realized through 
MA(^'O) [ig. In all cases the spectral function remains 
entirely self-consistent with the associated self-energy. 

Our test case for a momentum-independent self-energy 
is the simplest possible in momentum space, namely 
the Holstein polaron |34j : momentum- independent cou- 
pling between a single dispersionless phonon mode and 
tight-binding electrons. In reality however, even for 
the Holstein model, the self-energy is weakly dependent 
on momentum, which can be seen at the MA*^^^ level 
of approximation |32j . We overcome this by choosing 
the momentum-independent self-energy from the MA^^^ 
level in order to see how well these methods work for a 
truly momentum-independent self energy. For strongly 
momentum-dependent self-energy study we will model 
coupling to a single optical mode where the phonons live 
on half-integer lattice sites in between the electron sites 
and modify the on-site energy of their neighbours. In 
2D this describes lattice vibrations in a Cu02-like plane, 
where the motion of the O ions is the most important 
vibrational degree of freedom; this has been the topic of 
many ARPES studies [l0| .[2]M3] . Throughout the paper 
we will refer to this as the breathing-mode model. 

We may write both these models in the following form 
in momentum-space: 

k Q k,Q ^ 

(2) 

The terms describe, in order, an electron with disper- 
sion e^^ — —2t X^iLi cos(fcia) in D dimensions, an opti- 
cal phonon with energy D, and momentum Q, and the 
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FIG. 2: (Color online), (a) A{k,uj) calculated for the momentum-independent Holstein self-energy in ID with MA*^' for 
Q = 50meV and A — 0.5; the quasiparticle dispersion e', perturbation theory (about k = 0) prediction e^"^*^*, and bare- 
band are also shown, (b) Quasiparticle (vl), perturbation theory {v^'^^^), and bare-band (wjj) velocities, as well as (c) 
corresponding inverse masses, l/»n^, , and 1/m^ according to the definitions Vk = dsk/dk and l/m^ = d^Ek/dk^. 

(d) Momentum-dependent quasiparticle renormalization as obtained from v\/v1, ml/m'l, as well as the inverse quasiparticle 
strength l/Z^,, where Z'^ — J'^A{k,ij)duj is the quasiparticle-only integrated spectral weight; in the inset, these quantities are 
compared near fc = to the renormalization factors Q/W and (1 + A), obtained from quasiparticle bandwidth W (defined as 
the energy difference between top and bottom of the quasiparticle band) and dimensionless coupling X = /2tQ as well as 
the perturbation theory prediction for mass and velocity renormalizations (shown with the same line style, but which can be 
distinguished by their proximity to the quasiparticle curves) . In subsequent figures ([S] |4j [5| the quantity plotted is the effective 
coupling, Aeff, which would be implied by these renormalizations in the Migdal limit, which simply amounts to subtracting 1 
from the renormalization. (e-h) demonstrate similar traces for a momentum-dependent self-energy from coupling to a single 
breathing-mode in ID, here A = {\g\^) /2tQ, is the average of the coupling across the Brillouin Zone. The noise is due to the 
finite simulation grid and subsequent lineshape fitting; slight variations in peak position are enhanced by taking the derivative 
numerically and therefore most visible in m^. 



on-site momentum-dependent electron-phonon coupling 
gQ [for N sites with periodic boundary conditions; cj^ 
(ck) and (6^) are the usual electron and phonon cre- 
ation (annihilation) operators]. For the Holstein case 
gQ = g is a constant, leading to a dimensionless coupling 
X = g^ /2Dtil, the ratio between lattice deformation en- 
ergy —g^/^l and free-electron ground state energy —2Dt. 
For the breathing-mode gq = —i\/2g^-^^siii^Qia/2), 
which has an average value of (l^qP) = g^ across the 
Brillouin zone, allowing us to keep the same dimension- 
less coupling 49 . For this paper we set a = h = I and 
t — 50 meV, such that the ID bandwidth is 200 meV and 
the Brillouin zone is 27rA wide. Also note that an ad- 
ditional constant 1 meV FWHM Lorentzian broadening 
is used, similar to an impurity scattering, to allow the 
numerical resolution of the sharpest features in A{k,uj). 

The spectral function calculated with MA^^' for 
the Holstein problem in ID with ft = 50meV and 



A = 0.1,0.5,1.0 is presented as a false color plot in 
Fig.[lja,b,c), along with the path of peak maxima mea- 
sured through MDCs [/c„(w)] and EDCs (e^) compared 
with the bare-band dispersion e^. In the bottom panels 
(d,e,f) of Fig. [l] we present the real, S'(w), and imagi- 
nary, I]"(w), parts of the self-energy for this momentum- 
independent model. Here each is a true (and the low- 
est) pole of the Green's function (it has zero width, hence 
an infinite lifetime), and is only resolved numerically ow- 
ing to the impurity scattering inserted in the energy di- 
rection. One can see from Eq.[l] that the peak width 
should go roughly like E"(fc,w), and it is reassuring to 
see in Fig.[IJd,e,f) that the imaginary part of the self- 
energy is indeed zero near e^. The pole structure ef. is 
also distinct from that of km, the path of peak maxima 
observed during MDC analysis; not only are they funda- 
mentally different (as one is a function of uj and the other 
of k), but the path of peak maxima observed when cut- 
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FIG. 3: (Color online). Effective coupling, Aefi (as would 
be interpreted in the Migdal limit from the renormalization 
parameters defined as in Fig.[2|, plotted vs. the true dimen- 
sionless coupling A = {\g\^) /2tQ,; where g is a constant for 
the Holstein model (a), and qq = — i\/2<?sin(Qa/2) for the 
breathing-mode model (b). Also shown, in the inset only, are 
the predictions for observed effective coupling found via the 
mass renormalizations in perturbation theory (Eqs.js] and [4| 
for the low-coupling regime at A; = 0. Note that the noise in v 
and 1/Z a,t k — ki originates from the numerical determination 
of the inflection point ki. 



ting A{k, uj) in energy vs. momentum will not necessarily 
overlap, as has already been noted [il [55]. 

For small couplings [Fig. [T] panel (a)] most of the spec- 
tral weight remains along e^, with only a small feature 
formed at energy below the top of the band. With ex- 
perimental resolution such a feature might appear only 
as a "kink" in a quasiparticle dispersion, however from 
looking at the self-energy [panel (d)] one can see that a 
distinction between and km{^) should be made at this 
feature. The lowest pole, where S"(aj) ss and which we 
will identify as the quasiparticle, only exists between the 
top of the band and f2. This pole forms a narrow dis- 
persion, e^, of bandwidth approximately il, although for 
k near the zone edge the electron spectral weight is very 
weak due to it having significant phonon character. The 
km{i^) path of MDC peak maxima, however, does not 
follow this quasiparticle dispersion but instead carries on 
close to the original bare-band £\ into what we will iden- 
tify as the continuum due to its broader structure and 



finite E". 

As the coupling is increased [Fig. [T] panels (b,c)], this 
distinction becomes increasingly more evident; the quasi- 
particle band gains spectral weight toward the zone 
boundary and becomes more well defined. Also its band- 
width narrows, becoming less than $1 as the quasipar- 
ticle mass increases and the quasiparticle velocity de- 
creases. At the same time the spectral weight in the 
continuum becomes more spread out at deep energies, 
and new quasiparticle-like features begin to appear at 
the top of the continuum. At very large coupling (not 
shown) these additional features and the quasiparticle 
will eventually form a ladder of states with flat disper- 
sions, although this coupling regime is well beyond the 
scope of this paper. 



QUASIPARTICLE ANALYSIS 

As can be visualized from Fig.[l] quasiparticle renor- 
malizations do increase as the microscopic coupling in- 
creases. This monotonicity has led to widespread ac- 
ceptance of measuring coupling through the quasiparti- 
cle mass, velocity or strength renormalizations observed 
with ARPES, often without reference as to whether 
or not the system should be expected to fall in the 
Migdal/Eliashberg framework. In this section we use 
our simple models to demonstrate that this scheme 
is not universal, and to make other observations, by 
performing quasiparticle analysis as is typically done 
with ARPES data (Fig.[2]) on - 1200 generated spec- 
tral functions. These allow us to explore a wide range 
of couphngs (Fig. [3]), parameters (Fig.[4|, and different 
dimensionality (Fig.jsJ) on models which provide both 
momentum-dependent and momentum-independent self- 
energies. Following a discussion of these results we will 
follow the mass renormalization behavior as A — )■ for 
k ~ in detail through perturbation theory (see Fig. [6] 
and Eqs.|3]and|4|, the predictions for which are also plot- 
ted in Figs.[2]and[3]for comparison. 

In order to perform quasiparticle analysis we gen- 
erate an entire spectral function for each combina- 
tion of the parameters: model, dimensionless coupling 
A= (|(7|^)/2Z?i O, phonon energy scale ri/2t, dimension- 
ality D, and (for 2 and 3D) the desired cut through 
momentum space. For all simulations the form of the 
bare band is not changed and the hopping is set to a 
constant of i = 50meV to give physically familiar val- 
ues, a bandwidth of 200meV in the ID case. (To con- 
sider other bandwidths one should simply scale the band- 
width, phonon energy, and coupling together as seen in 
the Hamiltonian, Eq.[2|. On each of the ^ 1200 generated 
spectral functions the quasiparticle dispersion is found by 
fitting a Lorentzian peak with linear background to each 
EDC within the quasiparticle regime. The inclusion of 
a linear background allows the exclusion of any spectral 



5 




0.0 0.2 0.4 0.6 O.S 0.2 0.4 0.6 O.S 0.2 0.4 0.6 O.S 1 



X=-g^l2ta X=g^/2tn X=g^!2tn 

FIG. 4: (Color online). EfTective coupling, AcS, (as would be interpreted in the Migdal limit from the renormalization parameters 
defined in Fig.[2| from the Holstein model (a-c) and breathing-mode model (d-f) plotted vs. the true dimensionless coupling 
X^g^/2tn for a range of phonon energies Q labelled in meV. In panels (b, c) those from multiple renormalization parameters 
(which lie directly on top of each other) are all plotted, whereas in (e, f) those from inverse quasiparticle weight do not fall on 
any other curve and are therefore omitted for clarity. The slope at A = in panels (b, e) is the quantity plotted in Fig.lm Note 
that the noise in panels (c, f) originates mostly from the numerical determination of the inflection point fci, while in (a, d) it 
stems mostly from variations fitting the quasiparticle peak location at high momentum where it is has less weight. 



weight from the continuum which bleeds in (a problem 
especially at low couplings and high dimensions). We il- 
lustrate this analysis in Fig.[2j where we present the spec- 
tral function A{k,uj) for a mid-range coupling A = 0.5 
and phonon energy ft = 50meV for both the Holstein 
and breathing-mode models as well as dispersions found 
from the Lorcntzian fits e^, perturbation theory predic- 
tion e]^'^'^*, and the bare electronic structure [panels 
(a,e)]. Also shown are the velocities [w^, v^'^''* and in 
panels (b,f)], inverse masses and l/m^ 

in panels (c,g)], the corresponding renormalization ra- 
tios v\/v1, ml/m'l and their perturbation theory predic- 
tions, along with the inverse quasiparticle strength l/Z'^ 
and bandwidth renormalization il/W [panels (d,h), see 
caption for definitions]. In Fig. [2] panels (a-d) present 
the results for the Holstein model (with a momentum- 
independent self-energy), while panels (e-h) refer to the 
breathing-mode coupling (with a momentum-dependent 
self-energy). 

Fig. [2] (d,h) show that the velocity, mass, and spec- 
tral weight renormalizations are all functions of momen- 
tum, which raises concerns if one would compare them to 
the bandwidth renormalization fl/W^, or an 'expected' 



renormalization factor of (1 + A) which are both con- 
stant. Although they do cross at certain values of k this is 
merely accidental, and none match at the top of the band 
- our 'Fermi surface'. More problematic is that the mass 
renormalization must necessarily contain a divergence 
if the inflection point of e^, is different from e^, where 
vanishes (emphasized by the horizontal dashed 
line) . Similarly, in the case of momentum-dependent cou- 
pling [panels (e-h)], it can be seen that the quasiparticle 
dispersion si is not even monotonic, causing another di- 
vergence when vanishes in the middle of the dispersion 
(this non-monotonic dispersion is a direct consequence of 
the structure of the polaronic cloud, which causes a larger 
second-nearest-neighbour hopping, and is discussed at 
length in Ref. l36p . Due to this momentum-dependence, 
any estimation of A drawn from v^/v'J,, m'/m^., or 
would depend heavily on the momentum chosen; and if 
either of v\/vf, or rn1/m^^ were used close to their diver- 
gences, the estimated value could be off by an unlimited 
amount. Even il/W^, although constant in k, does not 
match the value of (1 -I- A) for either the momentum- 
independent or momentum-dependent case. From Fig. [2] 
we draw the conclusion that none of the renormalization 
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FIG. 5: (Color online). Effective coupling, Acfj, (as would be 
interpreted in the Migdal limit from the mass and velocity 
renormalization parameters defined in Fig.[2]| from the Hol- 
stein model and plotted vs. the true dimensionless coupling 
X = g^/2Dt Q for different dimensionalities and high-symmetry 
cuts. Note that the noise in v and 1/Z at k — ki originates 
from the numerical determination of the inflection point ki. 

quantities gives a good direct estimate of the dimension- 
less coupling X= {\g-^\) /2tQ. Further we find that, with 
the exception of the quasiparticle strength and velocity 
renormalization in the Holstein model only (which we 
will return to), the renormalizations do not match even 
each other - even though the models were kept as similar 
and simple as possible. This indicates that making even 
qualitative comparisons of 'coupling' from experiments 
on different materials (or even different experiments on 
the same material) through these renormalization param- 
eters may not be meaningful. However, modelling of the 
parameters in question from the original Hamiltonian via 
perturbation theory might be a start, as these results 
show much closer agreement near fc = despite the rel- 
atively high (for perturbation theory) coupling (we will 
return to discuss perturbation theory later). 

Despite their differences from each other and their mo- 
mentum dependence, however, these renormalizations do 
monotonically increase as the microscopic coupling in- 
creases (as previously observed in Fig.[T]) which naturally 
leads one to wonder how, precisely, these quantities scale 
with A as well as fl in our different models, so that one 
might be able to capture the trend if the material and 
measured quantity is held constant - for example in an 
experiment performed as a function of doping, if that 
doping does not cause structural distortions. In Fig. [3] 
we will follow the 'effective coupling', AcS (which is sim- 
ply the observed renormalization minus 1), that each of 
these renormalization quantities would predict using the 
Migdal/Eliashberg framework as a function of A, as well 
as renormalizations found using the perturbation theory 
results around k — (Eqs.[3]and|4|, for both momentum- 
independent [panel (a)] and momentum-dependent [panel 
(b)] self-energies. In Fig.|4] we plot a selection of these 
quantities in the same fashion, for a range of phonon en- 



ergies. For the momentum-dependent quantities we must 
choose a k value: we plot Vq/vq, mo/wQ, and 1/Zq at 
fc (our 'Fermi surface'), as well as v^./vf,. and 
at the inflection point k — ki of the quasiparticle band , 
where m'j.Jm\, diverges. 

In Fig. [3] we find that the predictions from all quan- 
tities scale monotonically with the microscopic coupling 
and are concave up. In the low coupling regime (below 
about A = 0.3) the perturbation theory results match the 
simulations - however nowhere does AeS match A from the 
Hamiltonian. At small coupling values, using this model, 
many renormalization quantities would drastically un- 
derestimate the true microscopic coupling, by a factor 
ranging from infinite {VL/W near A = where it is not 
renormalized in the Holstein model) to ~ 6 {\/Zk- near 
A = 0.1, breathing-mode). Conversely, at larger cou- 
pling values (A'-^2) all quantities would oiierestimate the 
true microscopic coupling, with factors ranging from ~4 
{^/W, Holstein) to ~22 (w^./w^., breathing- mode). We 
also find that, depending on the coupling or model, the 
relative renormalization strength of quantities changes 
- for the momentum-independent model VL/W is renor- 
malized the least, whereas in the momentum-dependent 
model shows the least renormalization. This indi- 

cates, yet again, that comparing different materials via 
renormalizations is not feasible, nor is comparing differ- 
ent renormalizations on the same material without a de- 
tailed model. We note again that, in the Holstein model 
only, quasiparticle strength and velocity renormalization 
are identical for all couplings at both fc = and k — ki 
(as previously seen in Fig. [2] where they are identical at 
all momenta). 

In Fig.g] we follow the same quantities for a variety of 
phonon energies, allowing Q, to vary from 1/16 to 3/4 of 
the bare-band width for both models (although inverse 
quasiparticle weight is omitted from the breathing-mode 
plots for clarity) . First we note that there are some qual- 
itative similarities, but just as many differences. In all 
these ID cases the concavity increases as phonon energy 
decreases, so that by the mid-coupling regime (A « 1) 
we recover the expected dependence - phonons which are 
easier to excite (require less energy) renormalize the band 
more. However, in the low coupling regime we do not find 
this dependence (later seen again in Fig. [6] and in agree- 
ment with Eqs.[3] and|4]). For both models the bandwidth 
[panels (a, d)] shows the opposite behavior for low cou- 
pling, with a transition near A = 0.8. Still considering 
the low-coupling regime, mass and velocity renormaliza- 
tions show little dependence on the phonon energy for 
the breathing-mode model, yet strong dependence in the 
Holstein case. Again we find that the renormalizations 
and their corresponding 'effective couplings' vary widely 
from each other, and depend on the model and parame- 
ters chosen - sometimes in counterintuitive ways. 

The final parameter to be varied is dimensionality, 
which we explore briefly with Fig. [5] for the Holstein 
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Hamiltonian in the low coupling regime only. Here we 
find that for a fixed dimensionality and phonon energy 
where ~ t the renormalizations as a function of A look 
qualitatively similar. The various renormalized quan- 
tities increase monotonically yet remain distinct from 
the microscopic coupling as well as each other (with the 
exception of quasiparticle strength and velocity renor- 
malization which are again the same), with details that 
depend on phonon energy and dimensionality. We feel 
it is important to note, however, that at larger cou- 
plings not explored here other studies on the dynamics 
of the Holstein (and momentum dependent Su-Schrieffer- 
Heeger) models have found more complicated behavior in 
higher dimensions; where a critical coupling value marks 
a drastic change in quasiparticle properties, which is most 
prominent as J7 — >■ [57^140) . However interesting, this 
type of behavior would not simplify quasiparticle renor- 
malization analysis on such a system and is not investi- 
gated here. 

So far we have shown that, while the slope may not be 
1, the renormalization curves could all still be reasonably 
well approximated as linear in A in the very low coupling 
regime and that perturbation theory makes an excellent 
prediction for them near fc = 0. This allows us to follow 
this slope more continuously through parameter space 
with perturbation theory than by simulating even larger 
numbers of spectral functions. It is worth noting that 
in the classic implementation of the (1 -I- A) scheme (see 
Ref.Jj, perturbation theory is discussed but dismissed as 
a possible avenue due to the resulting corrections being 
too large for perturbation theory to be valid. However, 
in that instance, some approximations are made to ease 
completion of the integrals which eliminate the possibility 
of the arbitrarily small couplings we have used here. In 
our case the lowest non-zero order in the phonon-electron 
interaction term from Eq.[2]is the second, and in ID we 
find that for both models it is possible to complete the 
integrals without further approximation. In the Holstein 
case the energy dispersion should be modified from £^ = 
— 2tcos(fca) to: 

£f * ^ -2t I cos(fca) + A ^ ) , (3) 

V ^(cos(fca) + g)'-i; 

which agrees with the results calculated for the quasipar- 
ticle residue at fc = in Ref.^T] For the breathing-mode 
model we find that the dispersion becomes: 



ITolstein Model (a) 
X Pertubalioii Theory 
• MA*'* Simulations 




Phonon energy (Qi'lt) 

FIG. 6: (Color online). Rate of change in observed effective 
coupling, Aoff, per change in true dimensionless coupling. A, 
defined as dXef[/dX\^_^ (where Acff = m^/m'l — 1 is as 
predicted by perturbation theory and measured from simu- 
lated spectral functions and A = (\g\^) /2DtQ,) plotted as a 
function of the phonon energy for (a) Holstein model and (b) 
coupling to a breathing mode. If the Migdal limit holds this 
derivative would be a constant, 1, for all values of Q/2t. For 
D > 1 cuts from the F point to both the corner and face of the 
Brillouin Zone were considered, as labelled. Vertical error in 
the results from simulated spectral functions is comparable to 
symbol size and results from approximating the slope at A = 
from the finite data points in Figs. |4] and [Sj as well as similar 
simulation sets not shown. Uncertainty in the perturbation 
theory results stems from the numeric Monte Carlo integra- 
tion used to determine the perturbation energies, taken from 
the distribution of independent subsets of total points evalu- 
ated, and is higher for non-diagonal cuts due to the narrower 
bandwidth in that direction. The curves predicted by Ref. l37l 
for the Holstein Model in 1 and 2 dimensions are not shown 
as they would be hidden by symbols, but fall exactly onto the 
perturbation theory results. 



« -2i(^cos(fca) + A|-F(i,fc)) 

/ON , sin^(fca) — ^cos(fca) , ^ 

J'{%,k) EE cos(fca) + — ^ ' ^* ^ ^. (4) 

y (cos(fca) + - 1 

This demonstrates that, at the very least, we should 
not expect the renormalizations to depend solely on the 



dimensionless coupling A — g^/2tn, but also on the other 
relevant energy scale in the problem - the ratio of phonon 
energy to bandwidth. By taking derivatives of these dis- 
persions we can also find the predicted mass and velocity 
renormalizations. In Fig. [2] we plot the predicted disper- 
sion, derivatives, and renormalizations contrasted against 
the simulated spectral function and find close but not 
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perfect agreement for both models near k = (but fail- 
ing at large momenta), despite the relatively high cou- 
pling (A — 0.5). As seen in Figs. |3] and |6] near A: = for 
vanishing A there is perfect agreement within our mea- 
surement accuracy; perturbation theory begins to show 
signs of failure near A « 0.25. In higher dimensions we 
did not complete the integrals exactly but instead used 
the VEGAS Monte Carlo integration algorithm to eval- 
uate them numerically Using this routine for 
all dimensionalities allowed us to validate the results by 
comparing them to these known solutions for ID for both 
models and the results found in Refs.[33 Elfor the Hol- 
stein model in 1 and 2D, where they show perfect agree- 
ment (once corrected for a slightly different definition of 
A in 2D). 

In Fig.|6] we show how the renormalization with the 
dimensionless coupling parameter A near fc = 0, A = 
(and hence the observed effective coupling, Acft) varies as 
a function of the other energy scale i7/2t for both models 
in 1, 2 and 3 dimensions, and how this matches nearly 
perfectly against measurements of the same quantity on 
the simulations. Interestingly, despite spectral functions 
which have outwardly similar characteristics [as seen, for 
example, in Fig. [2] panels (a vs. e) or Fig. 10 



Fig.[8], 



we find a drastic difference in how the actual renormaliza- 
tions vary with phonon energy depending on the model, 
and that neither model would be well approximated by a 
Acff = A scheme, which is shown as the dashed line along 
1. The ID Holstein model shows a perhaps expected de- 
pendence: phonons which require very little energy to 
excite have a dramatic effect on the electronic renormal- 
ization (blowing up as f2 — >■ 0); but as the phonon en- 
ergy increases, the mode has progressively less effect. In 
the 2D case, however, we find very limited dependence 
on phonon energy with a curve that is almost flat and 
could therefore be rescaled to match if A was chosen to 
be defined appropriately. In 3D we find the opposite 
of the ID case whereby the renormalization vanishes as 
— ?> 0. These very different limits are directly related to 
the bare-electron DOS at the band-edge, and its strong 
dependence on dimensionality |37j . In all dimensions we 
find that the renormalization is isotropic (as one might 
expect from an isotropic coupling) and that it asymp- 
totically approaches a similar value for large phonon en- 
ergies - reminiscent of a renormalization which depends 
solely on A, if only for r2/2t ^ 1. In stark contrast, 
however, the more realistic breathing-mode model shows 
entirely different behavior. For all dimensionalities the 
overall functional form is similar for — ^ 0, where mass 
renormalizations vanish. This low energy behavior may 
be expected as for weak coupling and energies close to 
A: = in the bare-band, the electron couples mostly to 
q ~ phonons; and in this model such coupling vanishes, 
g{q — !■ 0) — > 0. As the phonon energy increases, however, 
we discover that the renormalization is anisotropic with 
stronger renormalization along the diagonal cut (as may 



be expected for an anisotropic coupling) and a coupling 
which gets stronger as dimensionality increases (opposite 
the Holstein case) . We also find that the renormalizations 
do not asymptotically approach any fixed value for large 
phonon energies, as they did for the Holstein case. 

Overall we find that there is much variability in quasi- 
particle analysis, to the point that one cannot make a 
general rule about renormalizations in this regime. There 
are, however, two common threads. Firstly, for both our 
models, perturbation theory works in the low coupling 
regime: it correctly predicts the quasiparticle band struc- 
ture near fc = for all combinations of parameters tested, 
although it fails at higher momenta (as seen in Fig.[2|. 
The second, and perhaps more interesting, hints at some- 
thing which may be quantitatively gained through quasi- 
particle analysis - without even a more detailed model on 
which to attempt perturbation theory or other tools. In 
Figs.[2j |3j |4j and [5] one observes that, in the momentum- 
independent case only, v^/vl and ^/Zl lie precisely on 
top of each other for all values, and match rnj^/m\. at 
fc = 0. Although the velocity and mass renormalizations 
at fc = are simply a consequence of derivatives following 
each other near an extremum, the velocity renormaliza- 
tion and quasiparticle strength have implications for the 
structure of the self-energy, as was previously noted in 
Ref.i20i and is discussed in greater detail now. 

By definition, the Green's function is: 

G{k,u) = , , ^ . . (5) 

a; — — L(fc,a;) -I- ir) 

In the infinite lifetime quasiparticle regime the self- 
energy should be purely real, with any broadening com- 
ing solely from the small impurity scattering, ry. We may 
then identify the implicitly defined quasiparticle disper- 
sion as e| = e^-|-I](fc, w) and, since we are interested 

in an EDC, expand the self-energy about uj = e\ to first- 
order in energy. Taking — I/tt times the imaginary part 
yields the spectral function: 

^(^'^)^^^^ n..,^ 3S(fc,^) | 



■n 



2 + (c.-4)2(l- 



We can see that, cut along energy in the quasiparti- 
cle regime, the spectral function will be a Lorentzian 



at 
Zk 



with 
1/(1 



width 



given by 77 and with weight 
-I _ 5 ) . This relationship be- 



tween quasiparticle weight and the energy derivative 
of the self-energy is often derived, and usually associ- 
ated directly with the success of an effective coupling 
scheme [HI Ej, but we do not make that association 
here. Velocity renormalization is simply the ratio of mo- 
mentum derivatives of the bare, v' 



1 i 



dk 



f., and quasiparticle 
bands, which reduces 



to: 



1 



dk vl 



(7) 



We see that, for momentum-independent self-energies, 
the velocity renormalization must follow the inverse spec- 
tral weight. This means that the renormalization quan- 
tities can be used to conclusively check the momentum- 
dependence of the self-energy, in the quasiparticle regime. 
Whether or not the self-energy is momentum-dependent 
is of great importance to MDC self-energy analysis, on 
which we focus in the rest of the paper. 



SELF-ENERGY ANALYSIS 

Since quasiparticle analysis fails to reveal the true mi- 
croscopic coupling through renormalizations, we look to- 
ward other options for spectral function analysis. In this 
section we describe how it is possible to estimate the self- 
energy through the analysis of MDC lineshapes. We will 
begin, for simplicity, with a description for momentum- 
independent self-energy and move on to describe the 
implications of momentum-dependence on the proce- 
dure. Cases of a momentum-independent self-energy 
can be verified through quasiparticle analysis; as seen 
in Eq.j?] the overlap of v^/v'j, and '^/Zf, is only possi- 
ble when the self-energy is momentum-independent along 
the quasiparticle dispersion. Although the quasiparticle 
dispersion and the path of MDC peak maxima where 
MDC analysis is carried out may vary, in practice they 
are often very close in the quasiparticle regime. One 
must always keep in mind that although a momentum- 
independent self-energy causes a Lorentzian MDC line- 
shape, Lorentzian lineshape alone is not sufficient to con- 
clude that T.{k,uj) = T.{uj) [23\. 

In cases where the self-energy is momentum- 
independent we may analyze A(fc, a;) in terms of MDCs 
at constant energy uj — uj, where the self-energy may 
then also be considered a constant. Under this con- 
dition, as long as e\. can be linearized in the vicinity 
of the MDC peak maximum observed at = fcm, the 
MDC lineshape will be Lorentzian. By Taylor expand- 
ing e\. about an MDC peak maximum at fc = krm i-e. 
^\ ~ ^\ + ^fc ' {k — km) + ignoring higher order 
terms (which must be negligible if the curve appears 
Lorentzian) , and noticing that -|- — w = will 
implicitly define the observed peak maximum, we can 
rewrite Eq.[T]as: 



Aa,{ky. 



TT (A: - fc„,)2 + (Afc,„)2 



(8) 



with: 



Afc„, = =HWHM, 
Ao= l/vl = [Aa,{k)dk. 



(9) 
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Here Akm is the half-width half-maximum (HWHM) of 
a Lorentzian of weight ^o- For each constant energy cut. 



FIG. 7: (Color online). Eq.jsjand Eq. 10 described diagram- 
matically for a momentum independent self-energy, (a) is a 
false color plot with the bare-band (e^) and path of peak max- 
ima (km path) shown in addition to an example MDC cut at 
Cj and the peak maximum location for that cut [fcm(a})]. (b) 
the cut through momentum of A{k, uj) at constant energy Cj, 
observed to be a Lorentzian with a peak maximum located at 
fcm(a)), a HWHM Afcm(a)), and an area j4o(cD). These line- 
shape properties are related to the self-energy at uj through 
the bare-band evaluated at km through Eq.|10| 



Lu — UJ, the observed peak maximum is labelled km- The 
self-energies are then easily found as: 

j:'^ = -Akmvt- (10) 

One must be careful visualizing the relationships; al- 
though a simple picture might be that the peak, whose 
width scales with the imaginary self-energy and band ve- 
locity, has been 'pushed up' by the real self-energy to its 
observed location at w, one must remember that these 
quantities are defined implicitly and evaluated at differ- 
ent locations in the {k,uj) plane: the self-energy is eval- 
uated at {km,oj) and the bare-band at {km,£\ )■ Of 
course, these implicit definitions are no trouble if you 
simply wish to observe A{k, w) and not calculate it based 
on this simple approximation. These relations are illus- 
trated graphically in Fig.jTj 

These relationships work exactly where they are appli- 
cable: namely when the self-energy is momentum inde- 
pendent, km is far from a band edge (where v\ must 
vanish), where the peak shape is truly Lorentzian, and 
when the peak width is narrow enough that a first-order 
expansion of e\. is appropriate. Since the convolution of 
two Lorentzians is another Lorentzian where the peak 
width is a simple sum of the widths of the original func- 
tions, the inserted impurity scattering shows up directly 
as an addition to the measured (for comparison pur- 
poses a constant 77 = 1 meV has therefore been subtracted 
from all plots of E^). However, these relationships still 
hinge on knowledge of the bare-band. If e\ is unknown it 
is possible to fit it, to within an arbitrary energy offset, to 
any functional form which provides a value and derivative 
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FIG. 8: (Color online). (a) A{k,uu) calculated for the 
momentum-independent Holstein self-energy with MA'i' for 
n = 50meV and A = 0.1 as a false color plot; also shown are 
the km path of MDC maxima along which the analysis is per- 
formed, as well as the known bare-band and the third order 
polynomial approximation found through the KKBF analysis 
(the arbitrary energy offset introduced through KKBF has 
been shifted back by hand to allow comparison between the 
two). (b,c) Real and imaginary parts of the self-energy from 
the model {T,known), the bare-band and MDC fitting routine 

and the KK transform of Y^'Ijqq 
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(TiMDc) as found via Eq, 
{"^kk) and Ylj^jijQ {"^KK 
In (c) the MDC ratio results 
also shown. 



used as an internal check in KKBF. 
as found via Eq. 



T" 
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using a Kramers-Kronig bare-band fitting (KKBF) rou- 
tine (see Appendix). Alternatively, as previously noted 
in Ref. .20 and used in Ref. 45; the imaginary part of the 
self-energy requires knowledge of only u^.^ , which can be 
obtained directly from Aq — I /v\ , allowing us to write 
it as the MDC width/integral ratio: 

s;'„,„=-Afc„Mo (11) 



This variation allows us to tackle the problem over a 
larger range of A as it does not rely on the KKBF rou- 
tine to succeed over the entire range in order to ensure 
the fitness of the Kramers-Kronig transform and fit the 
bare-band. Eq.[TT]is free to work over energies where the 
peak is Lorentzian (i.e. Eq.|8]holds), and to fail in others 
without allowing these failures to block the procedure. 
Experimentally, when using this ratio, one must be care- 
ful that the observed spectra are properly normalized, 
otherwise it will be off by a constant multiple, but even if 
this is not possible the form of the imaginary self-energy 
should be nevertheless recoverable. It is also possible, 
in cases of momentum-independent self-energy for which 
v'l/v'l = 1/Z| from Eq.[7| to find the same ratio using only 

quasiparticle properties as I]'^£,c ~ "'^fc ^^m/Zl . 

The results of both the KKBF and the ratio method, 
performed as if the bare-band was not known on a 
momentum-independent self-energy in the low-coupling 
regime, are presented in Fig.jSj The internal self- 
consistency of the KKBF results is confirmed by the good 
match between Smdc and S^i^, and the agreement of 
'^Ratio adds further confidence. These 'measured' quanti- 
ties show good agreement with their known counterparts 
everywhere Eq.|8]is applicable, demonstrating that these 
methods work well in the low-coupling regime; however, 
they become progressively less accurate as the coupling 
increases. In Fig.[9]ja) we show the progressive failure of 
the method applied for couplings where A ranges from 
0-1, which demonstrates a rapid degeneration of the ac- 
curacy of the method outside of the low-coupling regime. 
Note, however, that the two methods fail in different 
ways. E^'j^p tends to fail more globally, whereas ^'^atio 
often continues to work almost exactly in some energy 
regions while failing catastrophically in others (these re- 
gions cause its average deviation, shown in Fig.[9j to in- 
dicate perhaps a higher degree of failure than deserved) . 
In Fig.[9]Jb) we demonstrate these differences by showing 
the results of the methods applied blindly at A = 0.5, 
even though lineshapes show that there are areas where 
Eq.[8]does not apply. One can see how the internal KKBF 
check has begun to fail as EX/dc and do not match; 
there are structural differences and E'^^ has picked up 
different offsets in the different flatter parts of the spec- 
trum, making setting its overall offset difficult. As the 
disagreement between E^'j^,^ and E^^ increases with 
coupling it will eventually cause the KKBF routine to 
fail to capture the bare electronic structure. None of the 
methods reproduce ^'pcnown "^^ar the sharp one-phonon 
structure; note that E"^^^ and ^'jiatio overestimate and 
underestimate it, respectively. Our experience with this 
model leads us to believe this to be typical behavior: 
when each method fails they do not tend to fail in iden- 
tical ways, so that in regions where they do match one 
can still have confidence that the methods are working. 

These methods hinge on the momentum-independence 
of the self-energy in two ways. For fitting line- 
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shape, a momentum-independent self-energy guarantees 
a Lorentzian lineshape but the inverse is not true - it is 
still possible to have a momentum-dependent self-energy 
which generates a Lorentzian. Additionally, fitting the 
bare-band relies on the Kramers-Kronig transforms in en- 
ergy, which are only valid for a fixed momentum. In cases 
where the momentum-dependence is not too strong lo- 
cally near km, however, we have found that it is still pos- 
sible to gain insight using similar approaches, although 
even more care must be taken in the interpretation of 
the results. If, despite momentum-dependence, the MDC 
appears Lorentzian in shape, certain higher order terms 
must be small when expanding each of e^, E'(a;,fc) and 
about km- Under this condition we may drop 
terms of order (k — k^)^ and higher from the denomina- 
tor as well as (fc — km) and higher from the numerator, 
which also implies we may drop dY/' /dk and d^TJ' jdk^ 
everywhere. We end up with an expression identical to 
Eq.|8j only with new interpretations for the HWHM as 
well as the spectral weight of the peak: 
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FIG. 9; (Color online), (a) Deviation (i.e., average of the 
squared difference at each km) between estimated and known 
self-energies vs. A. (b) Various estimates for the imaginary 
part of the self-energy, defined as in the caption of Fig.jS] for 
A{k,u) calculated within MA^^^ for J7 = 50meV and A = 0.5. 



^km 

Aci,{k)dk 



E"(lD,A:„ 



vl^^+dY'{Co,km)/dk' 
1 

vl +dY.'{Cj,km)ldk' 



(12) 



Now that the self-energies are momentum-dependent 
it becomes more important to remember that the self- 
energy extracted using this method will follow the km 
path through (w, k) space; this path is demonstrated as 
a false color plot in Fig. 10 |d,e,i,j). From Eq. 12 we find 
that, in this momentum-dependent case, the ratio check 
Eq.jTT] proves to be invaluable as it removes the possible 
inaccuracies when strong momentum-dependence near 
km in 5]' might provide a Lorentzian lineshape with a 
misleading width viewed from a momentum-independent 
perspective. In Fig. 10 'a-c,f-h) we present the results 
of both KKBF and the ratio method for a momentum- 
dependent coupling. From a comparison between the 
measured and known self-energies on paths through k 
space along the zone boundaries and along km, we find 
agreement only along the fc,„ path, as expected. Fur- 
ther, especially in Fig. |10| i near the 3rd phonon struc- 
ture close to ^ 125meV, one can see how it is pos- 
sible for Yj'lj^jj(j and S'^^ to agree with each other 
and yet not correctly predict 'S'Known^ despite the peak 
shape being reasonably Lorentzian, due to sufficient local 



momentum-dependence in the real self-energy (Fig. 10 



In this location we note, however, that S'^^^^^ still cor- 
rectly predicts ^Known it is not affected by this lo- 
cal momentum-dependence. Overall we find that, in 
a similar fashion to the momentum-independent case, 
there is generally good agreement between the found self- 
energies and the self-energy along the km path in the 
low-coupling regime and the methods progressively fail 



as we move into the mid-coupling regime. The available 
modes of failure are increased: there are more locations 
where the lineshape is not Lorentzian due to strong local 
momentum-dependence of the self-energy; places where 
it is still Lorentzian but with a misleading width; and the 
Kramers-Kronig relations are not valid along an arbitrary 
path through (w, fc) space, which disrupts the fitness of 
the KKBF routine. While it is not surprising that the 
KKBF routine eventually fails for large couplings in the 
momentum dependent case it is interesting that it works 
at all, as the Kramers-Kronig relations in energy are only 
formally valid for a fixed momentum but the measured 
self-energies follow the km path at all couplings. Despite 



this, as can be seen in Fig. 10 (b, c), the Kramers-Kronig 
relations appear to work relatively well along the km path 
in the low coupling case. Nevertheless we find that, in 
this model, failures occur at qualitatively similar cou- 
plings when momentum-dependence is added. 



CONCLUSIONS 

The spectral function is an extremely rich data source. 
We have shown that, despite its allure, away from the 
Migdal limit it is not always possible to extract the true 
microscopic coupling through quasiparticle renormaliza- 
tions of ARPES data with an effective coupling scheme 
- and certainly not for cases close to a dispersion max- 
imum. In this limit AeS ^ A. However, through care- 
ful modeling and the analysis of specific features it may 
be possible to learn much more. If one can correctly 
guess the model it may be possible to predict a given 
renormalization, or even show experimental support for 
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FIG. 10: (Color online). (a-c,f-li) A{k,Lj) and self-energies as defined in Fig.jSjfor momentum dependent coupling to a single 
breathing mode. In this instance, as the self-energy is momentum-dependent, the known self-energies are plotted along the 
path of peak maxima (SKnown@fcm) to show good agreement with the derived values, as well as along both edges of the 
Brillouin Zone for comparison (EKnown@fc — and EKnown@fc ~ tt). Finally (d,e,i,j) show the full extent of the energy and 
momentum-dependence of the real and imaginary self-energies as a false color plot, with the km path superimposed. 
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a given model via relationships between renormalization 
parameters. Through MDC self-energy analysis, we have 
shown that the self-energy can be extracted along paths 
through (fc,w) space in the low coupling regime - and 
it is possible still to gain some insight at higher cou- 
plings. We have also shown that it is possible to in- 
fer the momentum-dependence of the self-energy through 
comparison of renormalization properties. Methods like 
these, together with detailed models and powerful simu- 
lations, will hold the key to more thorough and quantita- 
tive analysis of the rich data supplied through ARPES. 



APPENDIX ON KKBF 

The method outlined here varies slightly from tech- 
niques previously described in the literature, which gen- 
erally reduce the functional form for substantially in 
order to facilitate an exact solution for A{k,u!) as they 
often deal with data very close to the Fermi energy over 
a small range [TH In our analysis we have in- 

stead expanded everything about km, essentially using a 
new linear approximation for on each MDC slice. Al- 
though our method has the disadvantage that it does not 
work as well near the zone boundary where the band ve- 
locity goes to zero (other methods which make a second 
order approximation can successfully predict and fit the 
non-Lorentzian shape in this region and may continue 
to work in this regime), ours has the distinct advantage 
that it works over a much larger energy range and allows 
fitting based on an infinite variety of bare-band models 
(so long as they are differentiable) . Most importantly, 
by its form it also explicitly shows that the self-energies 
are evaluated along the km path in the case where there 
is global momentum-dependence in the self-energy. One 
might imagine that for the analysis of a particular ex- 
periment one may have reasons to choose one method 
over another, or perhaps even a hybrid of the two. Here 
we will describe the idea of a Kramers Kronig bare-band 
fitting (KKBF) as implemented for our method, its ap- 
plication to other methods is similar. 

KKBF is a technique whereby a Lorentzian fit is first 
performed on every slice of constant energy, w, according 
to Eq.[8] The values of km and Akm from the fits can 
provide the self-energies for every {uj,km) point, within 
the limitations above, if the bare-band e^, is known. As an 
analytic complex function the real and imaginary parts 
of the self-energy are Kramers-Kronig related: 



Y,'^"^{k,u:)=±-V 

TT 



QC MDC 



(fc,C) 



(13) 



It is possible to 'fit' the bare-band parameters by 



choosing them such that E'^ 
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FIG. 11: (Color online). Self-energies as defined in Fig.js] 
for the MA'^' Holstein problem with Q, = 50meV and A = 
0.15, with extrapolated tails for Emdc and the KK transform 
shown. Panels (a) and (b) have the bias used in fitting the 
tails set too small, (e) and (f) have the bias set too large, and 
(c) and (d) have it set just right. 



^'kk- Since neither the Kramers Kronig relationships 



(Eq.|l3[) nor the MDC relationships (Eq. 10 1 are sensi- 
tive to a constant offset in both the real self-energy and 
bare-band this is unconstrained by the method and both 
S' and are free. In our study we have simply made 
the calculation of differences between "^kk and "Emdc 
insensitive to a constant offset, and set the final offset to 
zero by hand at the end to allow comparison. 

In our implementation of this idea, a simple third or- 
der polynomial was used to fit the bare-band with an ini- 
tial guess found by fitting MDC peak maxima. We then 
used the Levenberg-Marquardt Algorithm as imple- 
mented in the mpfit package for IDL [37] to vary band 
parameters. We found that the standard sum-of-squares 
minimization on the differences between Yikk and Em dc 
did not perform as well as a concave-down function, as 
it placed too much weight on outlying points far away. 
In order to evaluate the integrals in Eq. |13| with a finite 
region of data tails were extrapolated before a Fourier- 
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based transform was performed (the tails were then dis- 
carded, leaving the analysis of MDC and KK curves only 
within the data region). These tails were extrapolated 
by fitting an inverse polynomials to each side of the data, 
weighing the fit for each side's tail with an exponentially 
decaying bias parameter. A bias parameter of zero would 
weigh the entire curve equally, while a large bias param- 
eter would concentrate only near that data edge. 

It is possible for problems such as tails, overweighted 
outliers, and unconstrained offsets to compound each 
other. An unconstrained constant offset in E'^^^^ and 
e^. leads toward a tendency for a small linear offset in 
both, which when Kramers-Kronig transformed will dis- 
tort most visibly near the edges of the data, where 
it can interfere with a good fit of the tails. This in turn 
can lead to inaccuracies at these edges, which if over- 
weighted can distort the bare-band fit itself. This run- 
away condition results in a fit which gets progressively 
worse through iterations and will never find the correct 
bare-band. In practise we found that the tail bias param- 
eter as well as the the concavity of the function used to 
process errors must be carefully adjusted by hand in or- 
der to prevent this, which can be accomplished simply by 
looking at whether or not the tail approximation contin- 
ues to appear reasonable through successive iterations. 
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In Fig. 11 we show some typical examples of how the 
tail bias parameter can affect the fitting, each pair of 
panels represents the final "solution" of the entire band 
minimization problem using a given tail bias. Plots like 
these form the guide to be used when adjusting the bias 
parameters by hand while looking for the best solution. 
In panels (a) and (b) the tail bias is too small, and so 
the found tail is the best approximation which fits the 
whole curve. In panel (b) this causes a discontinuity for 
the low-energy tail right at the boundary, which in turn 
causes a cusp in the KK transform visible in (a). Despite 
this, the overall fit is not too bad with reasonable general 
agreement between MDC and KK self-energies - meaning 
the found bare-band is likely close to the real solution. In 
panels (c) and (d) the tail bias is good, which results with 
a realistic fit at all boundaries and a good agreement be- 
tween MDC and KK self-energies giving confidence that 
the found bare-band is accurate. In panels (e) and (f) the 
tail bias is too strong, which results in a tail fit depending 
too much on the data right at the edges. This results in 
a KK transform which is poor enough to thrown off the 
band fitting entirely resulting in a found bare-band which 
is likely not close to the true band, shown by generally 
poor agreement between MDC and KK self-energies. 
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